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Abstract — We consider a cyber-physical system consisting of 
two interacting networks, i.e., a cyber-network overlaying a 
physical-network. It is envisioned that these systems are more 
vulnerable to attacks since node failures in one network may 
result in (due to the interdependence) failures in the other 
network, causing a cascade of failures that would potentially 
lead to the collapse of the entire infrastructure. The robustness 
of interdependent systems against this sort of catastrophic failure 
hinges heavily on the allocation of the (interconnecting) links that 
connect nodes in one network to nodes in the other network. 
In this paper, we characterize the optimum inter-link allocation 
strategy against random attacks in the case where the topology of 
each individual network is unknown. In particular, we analyze the 
"regular" allocation strategy that allots exactly the same number 
of bi-directional inter-network links to all nodes in the system. 
We show, both analytically and experimentally, that this strategy 
yields better performance (from a network resilience perspective) 
compared to all possible strategies, including strategies using 
random allocation, unidirectional inter-links, etc. 

Keywords: Interdependent networks, Cascading failures, 
Robustness, Resource allocation, Random graph theory. 

I. Introduction 

Today's worldwide network infrastructure consists a web 
of interacting cyber-networks (e.g., the Internet) and physical 
systems (e.g., the power grid). There is a consensus that inte- 
grated cyber-physical systems will emerge as the underpinning 
technology for major industries in the 21st century J9|. The 
smart grid is one archetypal example of such systems where 
the power grid network and the communication network for its 
operational control are coupled together and depend on each 
other; i.e., they are interdependent. While interdependency 
allows building systems that are larger, smarter and more 
complex, it has been observed [ 17] that interdependent systems 
tend to be more fragile against failures, natural hazards and 
attacks. For example, in the event of an attack to an inter- 
dependent system, the failures in one of the networks can 
cause failures of the dependent nodes in the other network 
and vice versa. This process may continue in a recursive 
manner and hence lead to a cascade of failures causing a 
catastrophic impact on the overall cyber-physical system. In 
fact, the cascading effect of even a partial Internet blackout 
could disrupt major national infrastructure networks involving 
Internet services, power grids and financial markets [3]. Real- 



world examples include the 2003 blackout in the northeastern 
United States and southeastern Canada [17| and the electrical 
blackout that affected much of Italy on 28 September 2003 
0. 

A. Background and Related Work 

Despite recent studies of cascading failures in complex 
networks, the dynamics of such failures and the impact across 
multiple networks are not well understood. There is thus 
a need to develop a new network science for modeling 
and quantifying cascading failures, and to develop network 
management algorithms that improve network robustness and 
ensure overall network reliability against cascading failures. 
Most existing studies on failures in complex networks consider 
only the single network case. A notable exception is the 
very recent work of Buldyrev et al. [3] in which a "one-to- 
one correspondence" model for studying the ramifications of 
interdependence between two networks is set forth. This model 
considers two networks of the same size, say network A and 
network B, where each node in network A depends on one and 
only one node in network B and vice versa. In other words, 
each node in network A has one bi-directional inter-edge 
connecting it to a unique node in network B. Furthermore, 
it is assumed that a node in either network can function only 
if it has support from the other network; i.e., it is connected 
(via an inter-edge) to at least one functioning node from the 
other network. 

The robustness of the one-to-one correspondence model was 
studied in [3] using a similar approach to that of the works 
considering single networks Q, Q. Specifically, it is assumed 
that a random attack is launched upon network A, causing the 
failure of a fraction 1 — p of the nodes; this was modeled by a 
random removal of a fraction 1 — p of the nodes from network 
A. Due to the interdependency, these initial failures lead to 
node failures from network B, which in turn may cause further 
failures from network A thereby triggering an avalanche of 
cascading failures. To evaluate the robustness of the model, the 
size of the functioning parts of both networks are computed 
at each stage of the cascading failure until a steady state is 
reached; i.e., until the cascade of failure ends. One of the 
important findings of [3 1 was to show the existence of a critical 



threshold on p, denoted by p c , above which a considerable 
fraction of nodes in both networks remain functional at the 
steady state; on the other hand, if p < p c , both networks go 
into a complete fragmentation and the entire system collapses. 
Also, it is observed in [3 1 that interdependent network systems 
have a much larger p c compared to that of the individual 
constituent networks; this is compatible with the observation 
that interdependent networks are more vulnerable to failures 
and attacks. 

The original work of Buldyrev et al. |3j! has received much 
attention and spurred the study of interdependent networks in 
many different directions; e.g., see 0, 0, (TO], OH, lfT31l . 
|[T()1 . One major vein of work, including Q], lfT4l . Ifl6ll , aims 
to extend the findings of J3) to more realistic scenarios than 
the one-to-one correspondence model. More specifically, in [4 | 
the authors consider a one-to-one correspondence model with 
the difference that mutually dependent nodes are now assumed 
to have the same number of neighbors in their own networks; 
i.e., their intra-degrees are assumed to be the same. In lfT4ll the 
authors consider the case where only a fraction of the nodes in 
network A depend on the nodes in network B, and vice versa. 
In other words, some nodes in one network are assumed to be 
autonomous, meaning that they do not depend on nodes of the 
other network to function properly. Nevertheless, in lfT4l it was 
still assumed that a node can have at most one supporting node 
from the other network. More recently, Shao et al. |jT6 1 pointed 
out the fact that, in a realistic scenario, a node in network A 
may depend on more than one node in network B, and vice 
versa. In this case, a node will function as long as at least 
one of its supporting nodes is still functional. To address this 
case, Shao et al. |[T6l proposed a model where the inter-edges 
are unidirectional and each node supports (and is supported 
by) a random number of nodes from the other network. In a 
different line of work, Schneider et al. iflj | adopted a design 
point of view and explored ways to improve the robustness of 
the one-to-one correspondence model by letting some nodes be 
autonomous. More precisely, they assume that the topologies 
of networks A and B are known and propose a method, based 
on degree and centrality, for choosing the autonomous nodes 
properly in order to maximize the system robustness. 

B. Summary of Main Results 

In this study, we stand in the intersection of the two 
aforementioned lines of work. First, we consider a model 
where inter-edges are allocated regularly in the sense that all 
nodes have exactly the same number of bi-directional inter- 
edges, assuming that no topological information is available. 
This ensures a uniform support-dependency relationship where 
each node supports (and is supported by) the same number of 
nodes from the other network. We analyze this new model 
in terms of its robustness against random attacks via charac- 
terizing the steady state size of the functioning parts of each 
network as well as the critical fraction p c . In this regard, our 
work generalizes the studies on the one-to-one correspondence 
model and the model studied by Shao et al. [ 16|. From a design 
perspective, we show analytically that the proposed method of 



regular inter-edge allocation improves the robustness of the 
system over the random allocation strategy studied in [16|. 
Indeed, for a given expected value of inter-degree (the number 
of nodes it supports plus the number of nodes it depends upon) 
per node, we show that: i) it is better (in terms of robustness) to 
use bi-directional inter-links than unidirectional links, and ii) 
it is better (in terms of robustness) to deterministically allot 
each node exactly the same number of bi-directional inter- 
edges rather than allotting each node a random number of 
inter-edges. 

These results imply that if the topologies of network A 
and network B are unknown, then the optimum inter-link 
allocation strategy is to allot exactly the same number of 
bi-directional inter-edges to all nodes. Even if the statistical 
information regarding the networks is available; e.g., say it 
is known that network A is an Erdos-Renyi (2) network 
and network B is a scale-free network [1|, regular inter- 
edge allocation is still the best strategy in the absence of 
the detailed topological information; e.g., in the case where 
it is not possible to estimate the nodes that are likely to be 
more important in preserving the connectivity of the networks, 
say nodes with high betweenness |8|. Intuitively, this makes 
sense because without knowing which nodes play a key role 
in preserving the connectivity of the networks, it is best to 
treat all nodes "identically" and give them equal priority in 
inter-edge allocation. 

The theoretical results in this paper are also supported by 
extensive computer simulations. Numerical results are given 
for the case where both networks are Erdos-Renyi (ER) and 
the optimality of the regular allocation strategy is verified. 
To get a more concrete sense, assume that A and B are ER 
networks with N nodes and average degree 4. When inter- 
edges are allocated regularly so that each node has exactly 2 
bi-directional inter-edges, the critical threshold p c is equal to 
0.43. However, for the same networks A and B, if the number 
of inter-edges follows a Poisson distribution with mean 2, the 
critical p c turns out to be equal to 0.82. This is a significant 
difference in terms of robustness, since in the former case the 
system is resilient to the random failure of up to 57% of the 
nodes while in the latter case, the system is resilient to the 
random failure of up to only 18% of the nodes. 

To the best of our knowledge, this paper is the first work 
that characterizes the robustness of interdependent networks 
under regular allocation of bi-directional inter-edges. Also, it is 
the first work that determines analytically and experimentally 
the optimum inter-edge allocation strategy in the absence of 
topological information. We believe that our findings along 
this line shed light on the design of interdependent systems. 

C. Structure of the Paper 

The paper is organized as follows. In Section [II] we intro- 
duce the system model and present an overview of cascading 
failures. The behavior of this model under random attacks 
is analyzed in Section [Till where the functional size of the 
networks is characterized. Section[lV]is devoted to proving the 
optimality of regular inter-link allocation, while in Section [V] 



we give numerical examples and simulation results. Possible 
future research is explained and the paper is concluded in 
Section |VI] 

II. System Model 

We consider a cyber-physical system consisting of two inter- 
acting networks, say network A and network B. For simplicity, 
both networks are assumed to have N nodes and the vertex 
sets in their respective graphical representations are denoted 
by {vi, . . . , vn} and {v[, . . . , v' N }. We refer to the edges 
connecting nodes within the same network as intra-edges and 
those connecting nodes from two different networks as inter- 
edges. Simply put, we assume that a node can function only if 
it is connected (via an inter-edge) to at least one functioning 
node in the other network Q; and we will elaborate further 
on this. Clearly, the interdependency between two networks is 
intimately related to the inter-edges connecting them. In this 
study, inter-edges are assumed to be bi-directional so that it 
is convenient to use an N x N interdependency matrix C to 
represent the bi-directional inter-edges between networks A 
and B. Specifically, for each n, m = 1, . . . , N, let 

1 if v n and v' m depend on each other 
otherwise 



(C) 

nm 



(l) 



We also assume that inter-edges are allocated regularly so 
that each node has exactly k inter-edges, where k is an integer 
satisfying k < N. Without loss of generality, this strategy 
can be implemented in the following manner: For each n = 
1, 2, . . . , N, let the interdependency matrix be given by 



(C) t 



where we define 



if m — n,n i 
otherwise, 

n + e 

n + l- N 



.,ne (k - 1) 



(2) 



if 
if 



n + l < N 
n + l > N, 



for each I = 0, 1, . . . , N — 1; see also Figure Q] 

We are interested in evaluating the network robustness in the 
case of random node failures (or equivalently random attacks). 
Specifically, in the dynamics of cascading failures, we assume 
that a node is functioning at Stage i if the following conditions 
are satisfied J3], fl6l : i) The node has at least one inter-edge 
with a node that was functioning at Stage i — 1; ii) The node 
belongs to the giant (i.e., the largest) component of the sub- 
network formed by the nodes (of its own network) that satisfy 
condition i). For both networks, a giant component consisting 
of functioning nodes will be referred to as a functioning giant 
component. 

We assume that the cascade of failures is triggered by the 
failure of a fraction 1 — p of the nodes in network A. We further 
assume that these (1 — p)N nodes are chosen (say by the 
attacker) uniformly at random amongst all nodes in network 
A. By the definitions given above, it can be seen that after the 
initial attack, only nodes in the functioning giant component 
of A can operate properly. As a result of that, in the next 
stage, some of the nodes in network B may end up losing all 
of their inter-connections and turn dysfunctional. In that case, 



Network B 

4 




Fig. 1. A sketch of the proposed system model, namely the regular allocation 
strategy of bi-directional inter-edges: Each node in A is connected to exactly k 
nodes in B, and vice versa. 



TABLE I 

Key notation in the analysis of cascading failures 



A u Bi 


the functioning giant components in A and B at stage i 


PAi, PBi 


the fractions corresponding to functioning giant components 
at stage i, \Ai\ = PAiN, \Bi \ = p B iN 


Ai, Bi 


the remaining nodes in A and B retaining at least one 
inter-edge at stage i. 



the nodes that can function properly in network B will only 
be those in the functioning giant component of B. But, this 
fragmentation of network B may now trigger further failures 
in network A due to nodes that lose all their B-connections. 
Continuing in this manner, the cascade of failures propagates 
alternately between A and B, eventually (i.e., in steady state) 
leading to either: 1) residual functioning giant components 
in both networks, or 2) complete failure of the entire system. 
For an illustrative example, see Figure [2] where a cascading 
failure is demonstrated for a pair of interdependent networks 
with N — 6 nodes, k = 2, and p = 2/3. 

III. Analysis of Cascading Failures under 
Regular Allocation of Inter-edges 

In this section, we analyze the dynamics of cascading 
failures in two interacting networks. A principal objective 
of this study is to quantify the effectiveness of the regular 
allocation strategy for network robustness, by means of: i) 
characterizing the size of the remaining giant components in 
networks A and B after the cascade has reached a steady state, 
and ii) finding the corresponding critical threshold p c . To these 
ends, we will use the technique of generating functions fl2l . 
lfl3ll to analyze the sizes of functioning giant components in 
the two networks at each stage. For convenience, the notation 
used in the calculations is summarized in Table U 




Initial set-up Stage 1 Stage 2 Stage 3 Steady state 



Fig. 2. An illustration of cascading failures in two interdependent networks. Network A with nodes {vi , V2 , . . . , V e } and network B with nodes {v^v^, . . . ,v' 6 } 
are interdependent with each node having exactly two bi-directional inter-edges. Initially, a random attack causes the failure of nodes vi and V2- In stage 1, v\ 
and V2 are removed from the system along with all the links (inter and intra) that are incident upon them. As a result, node 113 becomes disconnected from 
the functioning giant component of network A, and thus fails. These failures then cause the nodes v' 2 and v' A to fail as they lose all their supports; i.e., all the 
inter-edges that are incident upon them are removed. In stage 2, we see the effect of removing v' 2 and v' 3 from network B: nodes and v' e fail as they become 
disconnected from the functioning giant component. The failure of nodes and v' 6 then leads to the failure of node v@ in stage 3, since was being supported 
solely byv'i and tig . By removal of the node vc, , the failures stop and the system reaches steady state. 



A. Stage 1 : Random Failure of Nodes in Network A 

Following the failures of a fraction 1 — p of randomly 
selected nodes in network A, the remaining network A\ has 
size pN; since we eventually let N grow large, pN can be 
approximated as an integer. As in |3|, 11121 . |fl3l , |[T6l , we use 
the technique of generating functions to quantify the fraction 
of the functioning giant component A\ C A\. Specifically, 
let the function Pa{p) determine the fraction of the giant 
component in a random subgraph that occupies a fraction p 
of the nodes in network A (the exact calculation of Pa{p) 
will be elaborated later). It follows that the functioning giant 
component has size 

\A 1 \=pP A (p)N := PA1 N. (3) 

As shall become apparent soon, at the end of each stage it 
is necessary to determine not only the size of the functioning 
giant component, but also the specific inter-edge distribution 
over the functioning nodes; i.e., the numbers of functioning 
nodes having particular numbers of inter-edges. Indeed, this 
is what makes the analysis of the regular allocation model 
more complicated than the models considered in [3], [14], 
i l76l/ . Here, at the end of Stage 1, each node in A\ still has k 
inter-edges from network B since network B has not changed 
yet. 

B. Stage 2 : Impact of Random Node Failures in Network A 
on Network B 

As the functioning part of network A fragments from A to 
Ai (in Stage 1), some of the inter-edges that were supporting 
S-nodes would be removed. Observe that the probability of 
removal can be approximated by 1 — |Ai|/|A| = 1 ~ pai for 
each inter-edge. With this perspective, a B-node loses k — j of 
its inter-edges with probability (j)p^i(l— pAi) k ~^ ■ Moreover, 
it stops functioning with probability (1 — pAi) k due to losing 
all k of its inter-edges. As a result, with B2 denoting the set 
of nodes in B that retain at least one inter-edge, we have 

\B 2 \ = (l-(l-p A1 ) k )N = p' B2 N, (4) 



where p' B2 = 1 — (1 — pAi) k ■ Also, the distribution of inter- 
edges over the nodes in B2 is given by 

\B 2 \ j =^jp j A1 (l-p A1 ) k -^N, j = l,2...,k, (5) 

with \B2\j denoting the number of nodes in B2 that have j 
inter-edges. 

As in Stage 1, the size of the functioning giant component 
B2 C B2 can be predicted by 

\B 2 \=p , B2 P B (p' B2 )N=p B2 N, (6) 

where Pb{) is defined analogously to the definition of P A (-) 
given in Section ITlI-AI Obviously, each node in B 2 can survive 
as a functioning node in B2 with probability Pb(p' B2 )- Thus, 
for each j — 1,2, ... ,k, the number of nodes in B2 that have 
j inter-edges is given (in view of (O) by 

\B 2 \j = P B {P B 2) Q)p> A1 (l-pAi) k - 3 N. (7) 

C. Stage 3 : Further A-Nodes Failures due to B-Node Failures 

Due to the fragmentation of the functional part of network 
B from B2 to B2 (not B to B2), some of the nodes in Ai 
may now lose all their inter-edges and stop functioning. To 
compute the probability of this event, first observe that each 
inter-edge from B2 to A\ will be removed with an approximate 
probability of 1 - \B 2 \/\B 2 \ = 1 - Pb(p' B2 )- Hence the 
probability that a node in A\ will lose all of its inter-edges is 
given by (1 — PB{PB2)) k ■ ^ a l so follows that the size of the 
network A3 C A\ comprised of the nodes that did not lose all 
their inter-connections is given via 

\A 3 \=p A1 (l-(l-P B (p' B 2)) k )N. (8) 

In other words, in passing from A\ to A3, a fraction 1 — 
|A 3 l/l Ai\ = (1~ PB(p'B 2 )) k of the nodes have failed. As pre- 
viously, the next step is to compute the size of the functioning 
giant component A3 C A3. However, this a challenging task 
as noted in [3|. Instead, we view the joint effect of the node 



failures in Stage 1 and Stage 3 as equivalent (in terms of the 
size of the resulting functional giant component; i.e., \A 3 \) to 
the effect of an initial random attack that targets an appropriate 
fraction (to be determined later) of the nodes. Intuitively, the 
node failures in Ai at Stage 3 (i.e., the removal of a fraction 
(1 — PB(PB2)) k °f nodes from Ai) have the same effect as 
taking out the same portion from A\ 0. In other words, it 
is equivalent to the removal of a fraction p(l — Pb{p'b2 )) k 
of the nodes from A. Recalling also that a fraction 1 — p of 
the nodes in network A failed as a result of the initial attack 
at Stage 1, we find that the fragmentation of A to A3 can as 
well be modeled (with respect to the size of A3) by an initial 
attack targeting a fraction 

1 - p + P (1 - P B {PB2)) k = 1 - P (1 " (1 - PB{ P ' B2 )) k ) 



of the nodes. It is now a standard step to conclude that, with 

PA3 = P (l — (1 ~ PB(PB2)) k )> me s i ze °f tne functioning 
giant component A3 is given by 



\A 3 \ =p' A3 P A (p' A3 )N =p A3 N. 



(9) 



D. Stage 4 : Further Fragmentation of Network B 

Due to the network fragmentation from A3 to A3 in Stage 
3, each inter-edge supporting a £>2-node will be disconnected 
with probability that equals the proportion nodes in A3 that did 
not survive to A 3 ; i.e., 1 - |A 3 |/|A 3 | = 1 - P A {p' A3 ) / P A {p) 
by ^ and (|9). Consequently, a node in B2 with j inter-edges 
will stop functioning with probability (1 — P A (p' A3 ) / P A {p)) 3 . 
Recalling also the inter-edge distribution ((T), the fraction L of 
node failures in B2 is given by 



1 vib I f^ Pa{p '^ 
^ |B2 H ^MpT 



p b(p B2 ) 1 -PAi 



A3) 



Pa{p) 



-d-PAi) 



Pb( P ' B2 ) {(1-pPa(p' A3 )) - (l~PAi) k 



Since \B 4 \ = \B 2 \ - LN, it follows that 

\B*\ = Pb( P ' B2 ) (l - (1 - pP A {p' A3 )) k ) N. (10) 

In order to compute the size of the functioning giant 
component B4 c B 4 , we proceed as in Stage 3. Specifically, 
we view the joint effect of node removals in Stage 2 and Stage 
4 as equivalent to that of an initial random attack which targets 
an appropriate fraction of the nodes. To determine this fraction, 
first observe that the failures in Stage 3 have triggered further 
node failures in B2 resulting a fraction 

1 - \B 4 \/\B 2 \ = 1 - (1 - (1 - P PA(p' A 3)) k )/p'B2 (ID 

of the nodes' failure. Next, note that the effect of these failures 
on \B 4 \ is equivalent to that of taking out the same fraction of 



nodes from B 2 |3|. Moreover, it has the same effect as taking 
out a fraction p' B2 {l - (l - (1 - pP A (p' A3 )) k ) /p' B2 ] of the 
nodes in B. Now, recalling that a fraction 1 —p' B2 of nodes in 
B have failed in Stage 2, we conclude that the joint effect of 
cascading failures in Stage 2 and Stage 4 (on \B 4 \) is identical 
to that of an initial random attack which targets a fraction 



PB2 



Pb 2 



1 - 



i-(i-pp a (p' A3 )Y 

PB2 



= (1 - pp A (p' A3 )f 

of nodes. As previously, with p' B4 = 1 — (1 — pP A (p' A3 )) k we 
conclude that the size of the functioning giant component B 4 
is given by \B 4 \ = p' B4 P B (p' B4 )N = p B4 N. 

E. Cascading Dynamics of Node Failures 

As mentioned earlier, the main goal of this section is to 
characterize the size of the functional giant components in 
steady state. Indeed, along the lines outlined above, one can 
obtain the sizes of all functioning giant components A\ D 
A 3 D ... D A 2m+ i and B 2 D B 4 D . . . D B 2m for any 
integer m. However, it is easy to observe the pattern in the 
expressions obtained so far and conclude that with p' A1 = p 
the size of all giant components are given by the recursive 
relations 



PA l =P A i P A(PAi), 
PAi=ph-0--PB(PBi-l))' 



3,5,7. 



(12) 



and 



PB l =P Bi Pb(p Bi )i oar n*\ 

, n h p ssfe i = 2,4,6,.... (13) 

This recursive process stops at an "equilibrium point" where 
we have p' B2m _ 2 = p' B2m and p' A2m -i = Pa^u+i so *at 
neither network A nor network B fragments further. Setting 
x = p' A2m+1 and y = p' B2m , this yields the transcendental 
equations 

x=p(l-(l-P B (y)) k ) y = 1 - (1 - P P A (x)) k . (14) 

The analysis carried out up to this point is valid for all 
networks, irrespective of their intra-structures. In principle, 
for specific intra-structures of networks A and B (which 
determine the functions P A and P B , respectively), the system 
(TPfl i of equations can be solved for given p and k. The steady- 
state fractions of nodes in the giant components can then 
be computed by using the relations lim.j_j.oo PAi '■— Pa^ = 
xP A (x) and lim^oc p Bi := P Boo = yP B (y). Indeed, in 
Section [V] we consider a special case where both networks 
A and B are Erdos-Renyi (ER) graphs [2 | and give solutions 
of the system (fT4l for several values of p and k. 

IV. Optimality of Regular Allocation Strategy 

In this section, we show analytically that the regular al- 
location strategy always yields stronger robustness than other 
strategies and thus it is optimal in the absence of intra-topology 
information. In the following, we refer to the system that uses 



the regular allocation strategy as System 1. Specifically, we 
consider two networks A and B where each node is uniformly 
supported by k bi-directional inter-edges. For convenience, 
we denote the fractions in the recursive relations ([T2li-([T3ll 
as PAi{p\ k) and p' Bi (p; k), where 1 — p is the initial fraction 
of failed nodes in network A. Also, we let Pa 1 ^ (pi k) and 
Pbi (p; k) be the steady-state fractions of functional giant 
components of the two networks, respectively. Finally, we use 
p Cl (k) to denote the critical threshold associated with System 
1. 

In what follows, we first investigate the dynamics of cas- 
cading failures in the auxiliary System 2, where bi-directional 
inter-edges are distributed randomly amongst nodes. The anal- 
ysis is carried out under a generic inter-degree distribution so 
that all possible (bi-directional) inter-link allocation strategies 
are covered. By making use of the convexity property and 
Jensen's inequality, we show that for a fixed mean inter-degree, 
System 2 achieves the highest robustness against random 
attacks when its inter-degree distribution degenerates, i.e.; 
when all nodes have exactly the same number of inter-edges 
so that System 2 is equivalent to System 1. Therefore, we 
conclude that regular allocation yields the strongest robust- 
ness amongst all possible (bi-directional) inter-link allocation 
strategies. Next, we show that systems with bi-directional 
inter-edges can better combat the cascading failures compared 
to the systems with unidirectional inter-edges 1161 . Together, 
these results prove the optimality of the inter-link allocation 
strategy in System 1; i.e., regular allocation of bi-directional 
inter-edges. 

A. Analysis of Random Allocation Strategy 

We now introduce the auxiliary System 2. Consider two 
arbitrary networks A and B, each with N nodes, and a discrete 
probability distribution F : N — >• [0, 1], such that 

F(j)=a j , j — 0,1, ... , (15) 

with Ej°lo otj = l. 

To allocate the interdependency links, we first parti- 
tion each network randomly into subgraphs with sizes 
aoN, ct\N, CK2N, .... By doing so, we can obtain subgraphs 

{S Aao ,SA ai ,S Aa2 , ■ ■ ■} and {S Baf> , S Bai , S Ba2 , ■ ■ ■}, such 
that 

\SA at \ = \S Bai \=ajN, j = 0,l,.... 

Then, for each j = 0,1,..., assume that each node in 
the subgraphs SA a and S Ba is assigned j bi-directional 
inter-edges. This ensures that the inter-degree of each node 
is a random variable drawn from the distribution F; i.e., an 
arbitrary node will have j inter-edges with probability aj, for 

each j — 0, 1, It is worth noting that the inter-degrees of 

the nodes are not mutually independent since the total number 
of inter-edges is fixed at E = ^ ctjjN for both networks. 

We have a few more words on the possible implementa- 
tion of the above random allocation strategy. Observe that 

'For N large enough, each of these subgraph sizes can be well approxi- 
mated by an integer. 



each bi-directional edge can be treated equivalently as two 
unidirectional edges. In this way, there are a total of 2E 
unidirectional inter-edges in the system, where E edges are 
going outward from network A and the other E edges are 
going outward from network B. We randomly match each 
unidirectional edge going outward from A to a unique edge 
going outward from B and combine them into a single bi- 
directional edge. To this end, let the edges going outward 
from A and B be separately labeled as e = {ei, . . . ,Ce} 
and e' = {e[, . . . , e' E }, respectively. Next, use the Knuth 
shuffle algorithm ifTTl to obtain random permutations e = 
{e\, . . . ,&e} and e! — {e^, . . . ,e' E } of the vectors e and 
e', respectively. Finally, for each i = 1, ...,23, match the 
unidirectional inter-edges e,j and e'j to obtain E bi-directional 
inter-edges. 

We now analyze the dynamics of cascading failures in 
System 2 using an iterative approach similar to that in Section 
lim For brevity, we skip most of the details and give only an 
outline of the arguments that lead to the sizes of functional 
giant components. The main difference from the analysis of 
Section HII] is that the fractions of nodes in A and B retaining 
at least one inter-edge, i.e., the fractions A, and Bi, need to 
be calculated differently from (0 and ( fTOb due to the random 
inter-degree of each node. 

Owing to the fragmentation from Bi-\ to £?i_i, each inter- 
edge supporting A could be disconnected with probability 1 — 
|2?j_i|/|2?j_i|, triggering further failures in network A at step 
i. With this insight, the aggregate effect of the failures in B 
up to stage i can be treated equivalently (with respect to the 
size of Ai) as removing each inter-edge supporting A with 
probability 1 — itj. According to Section [Till u% can be derived 
as follows: 

(i-D/2 . . 

U ^= II \w\ = p B{Pm-i) ^ = 3,5,7..., (16) 
£=i \ B ™\ 

Similarly, the aggregate effect of node failures in A before 
step i can be viewed as equivalent to removing each inter- 
edge supporting B with probability 1 — w, (with respect to the 
size of Bi) such that 

* = 77r n'jxH^^-i) * = 2,4,6.... (17) 

\ A \ l=1 \Au+l\ 

In System 2, each node is supported by j inter-edges with 
probability aj. In view of this, at step i, a node in network 

A would retain at least one inter-edge with probability 1 — 

00 

~ u i) 3 ■ Recalling also that a fraction 1 — p of the 

j=o 

nodes had already failed before the onset of the cascading 
failure, the equivalent remaining fraction of network A at stage 
i is given by: 

00 

PM = P(l -Yl a ^ X ~ U i) J ) 
3=0 

OO 

= P (i-j2<*j(i-PB(p'Bi-i)n 

j=0 



Similarly, the equivalent remaining fraction of network B turns 
out to be 



Proof: Since Pa and P B are monotonically increasing func- 
tions [13]. a sufficient condition ensuring (fJTJ will hold is 



3=0 



3=0 

Hence, the fractional sizes of the giant components at each 
stage are given (with p' A1 — p) by 



PAi =P'Ai P A(p' A i), 

Pm = P(l - E aiO-- p B(PBi-i)) j ), 
j=o 

for i = 3,5,7 ... , and by 

PBi = p' Bi Pb(p' Bi ), 

P'm = 1 -f:*j( 1 -P p A(p Ai -i)) i , 

3=0 



(18) 



(19) 



for i = 2, 4, 6, We next show that System 1 is always more 

robust than System 2 against random attacks by comparing the 
recursive relations (fT2l-([T3]l and (fT8l-([T9]>. 

B. Regular Allocation versus Random Allocation 

We now compare Systems 1 and 2 in terms of their robust- 
ness against random attacks. For convenience, we use a vector 
a. = (ao,Q!i, . . .) to characterize the inter-degree distribution 



F, where F(J) 



a 



Next, we denote the fractions in 



the recursive relations (fT8l-(fT9]> as p A i(p; a),p' Ai {jp; a) and 
PBiip; ct),p' Bi (p; a). Also, we let (p; a) and P B ^ (p; a) 
be the respective steady-state fractions of the functional giant 
components in the two networks where 1 — p is the fraction 
of initially failed nodes in network A. In other words, we set 
lim^oo p A i(p; ac) := Pjp^faa) and lim^oo p Bi (p; a) := 
P B 2 (p; a). Finally, we denote the critical threshold associated 
with System 2 by p C2 (a). 

Assume that network A (respectively network B) of Sys- 
tems 1 and 2 have the same size N and the same intra-degree 
distribution such that the functions Pa (respectively P B ) are 
identical for both systems. The next result shows that if the two 
systems are "matched" through their mean inter-degrees, i.e., 

oo 

if k = 2 atjj, System 1 always yields stronger robustness 

3=0 

than System 2 against random node failures. 
Theorem 4.1: Under the condition 



we have 



and furthermore 



PA^(p;k) > P A 2 o (p;o), 
p B^(p;k) > P B 2 o (p;a); 

P Cl {k) <p C2 (oi). 



(20) 



(21) 



(22) 



PAi(P\ k ) >PAi(Pi a )i < = 3,5,7..., 
p'm{p;k)>p' Bi (p;a), i = 2,4,6..., 



(23) 



where p' Ai (p;k), p' Bl (p;k), and p' Ai (p;a), p' Bi (p;a) denote 
the fractions in the recursive relations (fT2l-(fT3]). and (fT8l-(fT9l, 
respectively. We establish (l23l by induction. First observe that 
p' A iip; k) — p' A iip; a) = p and the inequality (l23l is satisfied 
for i = 1. In view of (QS-O and (TTBTl-(fl9Tl. condition © 
for i = 2 will be satisfied if 



oo 

{l-pP A (p' A1 [p;k)) k <Yai {l-pP A {p' AX {p; a )))\ 



J=0 



or equivalently 



(1 - P P A (p)) k < J2 a i 0- " pPa(p)Y- (24) 

3=0 

Under ( l20b . the convexity of (1 — pP A (p)) x implies (l24l by 
Jensen's inequality. Hence, we get that p' B2 (p; k) > p' B2 (p; a) 
and the base step is completed. 

Suppose that the condition ( l23l is satisfied for each i — 
1,2,..., 2m — 1, 2m. We need to show that d23l holds also 
for i = 2m + 1 and i = 2m + 2. For i = 2m + 1, the first 
inequality will be satisfied if it holds that 



(1 - PB(p' mm (j>; k))) k < E «3 (1 - PB(p' B 2m(P; «)))' 

3=0 

By the induction hypothesis, we have P B (p' B 2 m (P' ^)) — 
p B{p' B2m {Pi «)) since p' B2m (p; k) > p' B2m {p; a). As a result, 
the above inequality is satisfied if 

oo 

(l-u) k <^ Qj (l-«) j (25) 

3=0 

with u = P B (p' B2m (p;a.)). As before, under d20i i. ( f25b is 
ensured by the convexity of (1 — u) x in view of Jensen's 
inequality. The condition p' A2m +i{P^ k ) > P A 2m+i(P'> a ) is 
now established. 

Now let i = 2m + 2. The desired condition p' B2m+2 {p', k) > 
p' B2m+2 {p', ct) will be established if 

(i - pPA(p' A2m+ i(p; k ))) h 

OO 

3=0 



or equivalently 



(i-t0*<5>i(i-«) 3 ' 

3=0 



(26) 



where we set v = pP^(p^ 2OT , a)). The last step fol- 
lows from the previously obtained fact that p' A2m +i(P' k ) — 
p' A2m+ i(p', Oi) ■ Once more, d2oT ) follows by the convexity of 



(1 — v) x and Jensen's inequality. This establishes the induction 
step and the desired conclusion d2TT > is obtained. 

We next prove the inequality p Cx (k) < p C2 (a) by way of 
contradiction. Assume towards a contradiction that p c , 2 (a) < 
Pa{k) and fix p such that p C2 (ot) < p < p Cl (k). Then, let 
a fraction 1 — p of the nodes randomly fail in network A of 
both systems. Since p is less than p Cl , the node failures will 
eventually lead to complete fragmentation of the two networks 
in System 1; i.e., we get P A i (p; k) = P B i (p; k) = 0. On the 
other hand, the fact that p is larger than the critical threshold 
p C2 ensures P42 (p; a) > and P B 2 (p; a) > by definition. 
This clearly contradicts (f2Tb and therefore it is always the case 
that p Cl (k) < p C2 (a) under (|20|). ■ 

We have now established that the regular allocation of 
bi-directional inter-edges always yields stronger robustness 
than any possible random allocation strategy that uses bi- 
directional links. In the following section, we show that using 
bidirectional inter-edges leads to a smaller critical threshold 
and better robustness than using unidirectional inter-edges. 

C. Bi-directional Inter-Edges versus Unidirectional Inter- 
Edges 

We now compare the robustness of System 2 with that of 
the model considered in [16], hereafter referred to as System 
3. As mentioned earlier, the model considered in lfl6l is based 
on the random allocation of unidirectional inter-edges and 
can be described as follows. As with System 2, consider two 
arbitrary networks A and B, each with N nodes, and a discrete 
probability distribution F : N — > [0, 1] such that dl3T > holds. 
Assume that each node is associated with a random number of 
supporting nodes from the other network, and that this random 
number is distributed according to F. In other words, for each 
j = 0, 1, . . ., a node has j inward inter-edges with probability 
ctj. The supporting node for each of these inward edges is 
selected randomly amongst all nodes of the other network 
ensuring that the number of outward inter-edges follows a 
binomial distribution for all nodes. 

System 3 was studied in [ 16] using similar methods to those 
of Section [Til] and Section IIV-AI This time, after an initial 
failure of a fraction 1 — p of the nodes in network A, the 
recursive relations for the fractions of giant components at 
each stage turns [ 16 1 out to be (with p' A1 = p) 

PAi =p' Ai P A (Pm), 

p'm=p(i - E "ia-pfei-iWiH-i))''), (27) 

3=0 

for i = 3,5,7 ... , and 

PBi =Pm p B(Pm), 
P's 



i-E^-PM-iMPAi-i)) , 
3=0 



(28) 



fori = 2,4,6,.... 

Next, we compare System 2 and System 3 using the 
recursive relations (fT8l-JT9l) and (|2~7Ti-(|28]i. In doing so, we 
use the same notation to define the fractions in the recursive 
relations (fT~8b-(fT9b as used in Section ITV-Bl while the fractions 



in (|2~7Ti-(|28]i will be denoted by p Ai (p; a) , p' Ai (p; a) and 
Psi(p;")»P'si(P; a )- We let P A ^{p;a) and P B ^(p;a) be 
the steady-state fractions of functional giant components in 
System 3 if a fraction 1 — p of the nodes initially fail in network 
A. In other words, we set lim^oo p Ai (p; a) := Pa^ (p;ct) 
and lim^oo p Bi (p; a) :— P B s (p;a). Finally, we denote by 
p C3 (a) the critical threshold for System 3. 

The next result shows that System 2 is always more robust 
than System 3 against random node failures. 

Theorem 4.2: We have that 



and furthermore, 



P A 2 e (p;o) > P A ^(p;a), 
P B 2 (p;a) > P B3 (p;a), 



Pc 2 (oc) <p C3 (a). 



(29) 



(30) 



Proof: Since Pa (x) and P B (x) are monotonically increasing 
IPT31 . a sufficient condition ensuring (|29b is given by 



P A i(p; a ) >p'Ai(p; a )> i = 1,3,5.. 
PBi(P'> a ) >P'm(p; a ), i = 2,4,6. 



(31) 



We establish OTb by induction. First, observe that for i = 1, 
p' A1 (p;a) = p' A1 (p;ct) = p and condition (|3TT > is satisfied. 
Next, for i = 2, we see from (fl~9b and (f28b that the inequality 

p B2 (p;at) >p' B2 (jp;a) 

will hold if 

00 

^^(l-pP^foa)))^ (32) 

j=o 

OO 

< J2 a 3( l -PAl(P; a ) P A(PAl(P; a ))y- 

0=0 

Since p' A1 (p;a) = p' 3 A1 (p;a) = p, it is immediate that (|32| > 
is satisfied with equality and this completes the base step of 
the induction. 

Suppose now that condition (l3Tb is satisfied for all i = 
1, 2, . . . , 2m - 1, 2m. We will establish (ED for i = 2m + 1 
and i = 2m + 2 as well. Comparing (fT8l and (|27| |. it is easy 
to check that for i = 2m + 1, (l3~TT l will hold if 



^aj(l - P B (PB2 m (p; a ))Y 

3=0 



(33) 



< 



OO 



£^ M 1 - P'B2rn(P; °*) P B {p' B2m (p\ a W ■ 

3=0 



By the induction hypothesis, ( BTT l holds for i = 2m so that 
p B(p' 3 B2 m (j>;ot)) < P B (p' B2m (p;a)). It is now immediate 
that (l33l holds, since we always have p' B2m (p;a) < 1. This 
establishes OTT l for i = 2m + 1; i.e., that 



P'A2m+l{P: a ) < PA2m+l(P; ") 



(34) 



For i = 2m + 2, we see from (fT9] l and d28l ) that condition 
OTb will be satisfied if 



- pP A (p' A 2m+l(P- a W 



(35) 



J=0 



A2m+1 A2m+1 

3=0 

In view of (l34l and the fact that p'^2m+i(Pi a ) — P> 
we immediately obtain ( TJBI l and the induction step is now 
completed. This establishes condition d3TT > for alH = 1, 2, . . ., 
and we get d29l ). 

The fact that d29b implies (f30b can be shown by contradic- 
tion, as in the proof of Theorem 14.11 | 



Summarizing, it can be seen from Theorem 14.21 that using 
bi-directional inter-edges (System 2) always yields stronger 
system robustness compared to using unidirectional inter- 
edges (System 3). This being valid under an arbitrary dis- 
tribution a of inter-edges, we conclude that regular allocation 
of bi-directional inter-edges leads to the strongest robustness 
(amongst all possible strategies) against random attacks as we 
recall Theorem 14. 11 



V. Numerical Results: The Erdos-Renyi Networks 

Case 

To get a more concrete sense of the above analysis results, 
we next look at some special cases of network models. In 
particular, we assume both networks are Erdos-Renyi networks 
0, with mean intra-degrees a and b, respectively. For this 
case, the functions P A (x) and Pb{v) that determine the size 
of the giant components can be obtained [13] from 

P A (x) = l-f A and P B (y) = l-f B , (36) 

where f A and f B are the unique solutions of 

f A = exp{ax(f A -l)} and f B = exp{by{f B - 1)}. (37) 

In what follows, we derive numerical results for the steady- 
state giant component sizes as well as critical p c values. 
Specifically, we first study System 1 by exploiting the recursive 
relations (fT2l-(fT3l using ( f36l > and ( 1371 ). Similarly, we derive 
numerical results for System 2 by using the recursive relations 
<TT~8T> - (TT~9b . For both cases, we use extensive simulations to 
verify the validity of the results obtained theoretically. 

A. Numerical Results for System 1 
Reporting d36l > into (fl4l) . we get 

x = p{l-f k B ) y = l-(l-p(l-f A )) k . (38) 

It follows that the giant component fractions at steady state 
are given by 



p Aoo =p(i-/£)(i-/a), 

Pn=(i-a- 



(l-(l-p(l-/ A )) fe ) 



(39) 




d) p=0.563 



f) p=0.7 



Fig. 3. Possible solutions of the system (41t are depicted for several different 
p values when a = b = 3 and k = 2. In all figures, the x-axis represents Ja 
while the y-axis represents fg. The critical p c corresponds to the case where 
there is only one non-trivial solution to the system, i.e., the case when the two 
curves are tangential to each other. 



Next, substituting (l38l l into (|37| > we obtain 

f A = exp{ap(l - - 1)}, 

f B = cxp{6 (l - (1 -p(l - f A )) k ) (f B - 1)}. <40) 

We note that the system of equations d40b always has a trivial 
solution f A = f B = 1, in which case the functional giant com- 
ponent has zero fraction for both networks. More interesting 
cases arise for large values of p when there exist non-trivial 
solutions to d40b . In particular, we focus on determining the 
critical threshold p c ; i.e., the minimum p that yields a non- 
trivial solution of the system. Exploring this further, we see 
by elementary algebra that ( f40b is equivalent to 



log fA 

(fA-l)ap 



if Q < f A < i; v/ B if f A = i 



fA = 1 - 



'°B fB 

' (f B - 1 ) b 



V 



if < /b < 1; V/ A if f B = !• 

(41) 

In general, it may be difficult to derive an explicit expression 
for p c . Instead, we can solve ( PiTT l graphically for a given set 
of parameters a, b, k,p and infer the critical threshold p c using 
numerical methods. For instance, Figure [3] shows the possible 
solutions of the system for several different p values when 
a = b = 3 and k — 2. In Figures 0a-c), we have p < p c 
and there is only the trivial solution f A = f B = 1 so that 
both networks go into a complete fragmentation at steady 
state. In Figure Od), we have p = p c and there exists one 
non-trivial solution, since the two curves intersect tangentially 
at one point. In Figures Oe-f), we have p > p c and there 
exist two non-trivial intersection points corresponding to two 
sets of giant component sizes. In these cases, the solution 
corresponding to the cascading failures should be the point 
that yields the larger giant component size. In other words, 
the solution corresponds to the intersection point that is closer 
to the starting point of the iterative process (see d39l)). 



In the manner outlined above, we can find the critical 
threshold p c for any fixed values of the parameters a, b and 
k. As illustrated in Figure [3] we can further add the tangential 
condition 

f^X§L = l (42) 
djB djA 

to the equations ( |4TT > since the critical p c value corresponds to 
the tangent point of the two curves given by ( HTt . Thus, the 
critical values f Ae , Jb c an d Pc can be computed (numerically) 
for any given set of parameters through the following system 
of equations: 



and 



Sb = \1-TT^£- X0<f A <l; 
V (f A - l)ap 



fA = 1-- 

dfA 
df B 



)ap 

log/g 

(fB-l)b 



(43) 



if < f B < 1; (44) 



Eq. 



z.z. 



dh 
dft 



Eq. 



33 



= l. 



(45) 



The analysis results are now corroborated by simulations. 
In Figure |4(a)| we show the variation of p c with respect to 
k for different values of a = b, where the critical p c values 
are obtained by solving the system (|43T > graphically. To verify 
these findings, we pick a few sets of values a, b and k from 
the curves in Figure |4(a)| and run simulations with N = 5000 
nodes to estimate the probability p m f of the existence of a 
functional giant component in steady state. As expected 0, 
in all curves we see a sharp increase in p m f as p approaches 
a critical threshold p c . It is clear that the estimated p c values 
from the sharp transitions in Figure |4(b)| are in good agreement 
with the analysis results given in Figure |4(a)| 

B. Numerical Results for System 2 

As in System 1, the recursive process (fT~8T> - dT9b of System 
2 stops at an "equilibrium point" where we have p' B2m _ 2 = 

P' B 2m = x and P'A2,n-l = PA2m+l = V- This y ields the 

transcendental equations 



x= P (i- £ aj (i - P B {y)) 3 ), 

OO 

3=0 



(46) 



The steady-state fraction of nodes in the giant components can 
be computed by using the relations lim^oo p A % ■— P Aoa — 
xP A {x) and Yiu\ i ^r X p Bl := P Baa = vPb{v)- 

In particular, we assume that the inter-degree distribution F 
at each node is a Poisson distribution with mean k, and hence 



k J 

e~ k — , j = 0,1, 2, ...,oo. 

r- 



Substituting (O and (|47j into (|46j>, we get 



x = p \ 1 — 



j=0 



e~ k fB 3 I = p ( 1 



-k{l-f B ) 



(47) 



(48) 



l-ElT^^-Kl-ZA))^! 



-k P (i-f A ) 



■ (49) 



3=0 



Next, putting (gHJl and (|49j into (O, we find 



/A = l + ^ln(l 
Sb = 1 + \ In (l ■ 
V/a if /b = 1 



i^feVtf 0</b<1; 

+ s^)>if o</x<i; 

V/b if Sa = 1. 



(50) 



As in the case for System 1, the critical threshold p c for System 
2 corresponds to the tangential point of the curves given by 
d50l l. and can be obtained by solving d50b graphically. 

We now check the validity of these analytical results via 
simulations. In Figure |5(a)| we show the variation of analyt- 
ically obtained p c values with respect to average inter-degree 
k for different values of a = b. To verify these results, we 
pick a few sets of values a, b and k from the curves in 
Figure |5(a)| and run simulations with N = 5000 nodes to 
estimate the probability p m f of the existence of a functional 
giant component in steady state. As expected J5], in all curves 
we see a sharp increase in p m f as p approaches a critical 
threshold p c . It is also clear from Figure |5(b)| that, for all 
parameter sets, such sharp transition occurs when p is close 
to the corresponding p c value given in Figure |5(a)| 

C. A Comparison of System Robustness 

In Section IIV-BI and IIV-CI we have analytically proved 
that the regular allocation of bi-directional inter-edges leads 
to the strongest robustness against random attacks. To get 
a more concrete sense, we now numerically compare the 
system robustness of these strategies in terms of their critical 
thresholds p c . Specifically, we consider coupled Erdos-Renyi 
networks with mean intra-degrees a and b. For the sake of 
fair comparison, we assume that the mean inter-degree is set 
to k for all systems; in both Systems 2 and 3, the inter-degree 
distribution F at each node is assumed to be Poisson. The 
critical threshold value p c corresponding to all three strategies 
are compared under a variety of conditions. For Systems 1 
and 2, we use the numerical results derived in Section IV-AI 
and Section IV-BI respectively, while for System 3 we use the 
numerical results provided in |16|. 

First, we compare System 1 with System 3 to see the 
difference between the proposed regular inter-edge allocation 
strategy and the strategy in [16|. Figure [6(a)| depicts p c as a 
function of mean inter-degree k for various values of a = b, 
while Figure |6(b)| depicts the variation of p c with respect to 
a = b for different k values. In all cases, it is seen that regular 
allocation of bi-directional inter-edges yields a much smaller 
p c (and thus, a more robust system) than random allocation of 
unidirectional inter-edges. For instance, for a = b = k = 4, 
System 3 |fTi5l Figure 2] gives p c = 0.43, whereas, as seen via 
Figure [6(a)l System 1 yields a critical threshold at 0.317. This 
is a significant difference since it means that System 3 can 
have a functioning giant component despite a random failure 
of at most 57% of the nodes, whereas System 1, which uses the 
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(a) (b) 

Fig. 4. a) The critical p c value versus k for the regular allocation strategy (System 1 ). The plots are obtained by solving the system f45D graphically for various 
a, b values. It can be seen that as k increases the robustness of the system increases and the critical fraction p c approaches that of a single network; i.e., — b) 
Experimental results for the regular allocation strategy (System 1) with N = 5000 nodes. A fraction 1 — p of the nodes are randomly removed (from network 
A ) and the corresponding empirical probability p m f for the existence of a functional giant component at steady state is plotted. As expected, in all cases there is 
a sharp increase when p approaches a critical threshold p c ; for (a = b = 3, k = 3), (a = b = 3, k = 5) and (a = b = 6, k = 3), the critical p c values are 
roughly equal to 0.47, 0.41 and 0.23, respectively. Clearly, these p c values are in close agreement with the corresponding ones of Figure \4(a)\ which are obtained 
analytically. 
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(a) (b) 

Fig. 5. a) The critical p c value versus k for the random allocation strategy (System 2). The plots are obtained by solving the system (50b graphically for various 
a, b values. It is seen that the critical p c can be larger than one in some cases (e.g., for a = b = 3 and k = 1) meaning that the system collapses already without 
any node being attacked. This is because, due to the random allocation of inter-edges, a non-negligible fraction of the nodes receive no inter-edges and become 
automatically non-functional even if they are not attacked, b) Experimental results for System 2 with N = 5000 nodes. A fraction 1—pof the nodes are randomly 
removed (from network A) and the corresponding empirical probability p m ffor the existence of a functional giant component at the steady state is plotted. As 
expected, in all cases there is a sharp increase when p approaches to a critical threshold p c ; for (a = b = A, k = 2), (a = b = 4, k = 3) and (a = b = k = A), 
the critical p c values are roughly equal to 0.480, 0.380 and 0.335, respectively. Clearly, these p c values are in close agreement with the corresponding ones of 
Figure \5(a)\ which are obtained analytically. 



regular inter-edge allocation scheme proposed in this paper, is 
resistant to a random failure of up to 68% of the nodes. Indeed, 
in some cases, our strategy can outperform that in lfl6l even 
with half the (mean) inter-degree per node. For instance, when 



a = b — 4, our strategy yields p c = 0.414 with only k = 2 as 
compared to p c = 0.43 of the System 3 with k = 4. 

We also compare System 1 with System 2 in order to see 
the improvement in allocating bi-directional edges regularly 



S - a, b=3— System 3 
• — a,b=3 — System 1 
K - a, b=6— System 3 
A — a,b=6 — System 1 




■ ▲ - k=2 — System 3 

■ K - k=3— System 3 

■ • - k=4— System 3 
-A — k=2 — System 1 

-X — k=3 — System 1 
k=4 — System 1 
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Fig. 6. A comparison of System 1, System 2 and System 3 in terms of their critical p c values when expected inter-degree of any node is set to k. For System 2 
and System 3, the distribution of the number of inter-edges is assumed to be Poisson. In all figures, dashed lines correspond to System 3, dash-dot lines represent 
System 2, and solid lines stand for System 1. a) p c v.s. k is depicted for different values of a = b in System 1 and System 3. b) p c v.s. a = b is depicted 
for various k values in System 1 and System 3. In all cases, we see that the regular allocation of bi-directional inter-edges yields a smaller p c than the Poisson 
distribution of unidirectional inter-edges with the same mean value k. c) p c v.s. k is depicted for different values of a = b in System 1 and System 2. It is clear 
that System 1 yields a lower p c (and thus a higher robustness) than System 2 in all cases, d) p c v.s. a = b is depicted for various k values in System 1, System 2 
and System 3. In all cases System 1 yields the lowest p c (i.e., highest robustness), while System 3 has the highest p c (i.e., lowest robustness) and System 2 stands 
in between. 



rather than randomly. Figure |6(c)| depicts p c as a function of 
mean inter-degree k for various values of a = b. It is seen 
that, in all cases, System 1 yields a lower p c (and thus a more 
robust system) than System 2. For example, when a = b = 3 
and k — 2, we get p c — 0.56 for System 1, while for System 
2, we find that p c — 0.68. The difference is significant in that 
it corresponds to a resiliency against a random failure of up to 
44% of the nodes in System 1 as compared to 32% in System 
2. 

Finally, in order to better illustrate the optimality of System 
1 in terms of system robustness, we depict in Figure |6(d)| 



the variation of p c with respect to a = b for different values 
of k in all three systems. It is clear that the proposed regular 
allocation strategy in System 1 always yields the lowest p c and 
thus provides the best resiliency against random attacks. We 
also see that System 2 always outperforms System 3, showing 
the superiority (in terms of robustness) of using bi-directional 
inter-edges rather than unidirectional edges. 

We believe that the drastic improvement in robustness 
against random attacks seen in System 1 has its roots as 
follows. First, in the absence of intra-topology information, 
it is difficult to tell which nodes play more important roles 



in preserving the connectivity of the networks. Thus, in order 
to combat random attacks, it is reasonable to treat all nodes 
equally and give them equal priority in inter-edge allocation. 
Secondly, in Systems 2 and 3, there may exist a non-negligible 
fraction of nodes with no inter-edge support from the other 
network. Those nodes are automatically non-functional even 
if they are not attacked. But, the regular allocation scheme 
promises a guaranteed level of support, in terms of inter- 
edges, for all nodes in both networks. Finally, using bi- 
directional inter-edges ensures that the amount of support 
provided is equal to the amount of support being received for 
each node. Thus, the use of bi-directional inter-edges increases 
the regularity of the support-dependency relationship relative 
to unidirectional inter-edges, and this may help improve the 
system robustness. 

VI. Conclusion and Future Work 

We study the robustness of a cyber-physical system in which 
a cyber-network overlays a physical-network. To improve 
network robustness against random node failures, we develop 
and study a regular allocation strategy that allots a fixed 
number of inter-network edges to each node. Our findings 
reveal that the proposed regular allocation strategy yields the 
optimal robustness amongst all strategies when no information 
regarding the intra-topologies of the individual networks is 
available. For future work, we conjecture that in the presence 
of such information, the topology of the networks can be 
exploited to further improve the robustness of cyber-physical 
systems against cascading failures. 

It is also of interest to study models that are more realistic 
than the existing ones. For instance, in a realistic setting, one 
can expect to see a certain correlation between the inter-edges 
and the intra-edges of a system owing to the geographical 
locations of the nodes. Also, some of the nodes may be 
autonomous, meaning that they do not depend on nodes of 
the other network to function properly; in that case, one can 
expect the regular allocation strategy to still be the optimum 
if the nodes that are autonomous are not known. Clearly, 
there are still many open questions centered around network 
interdependence in cyber-physical systems. We are currently 
investigating related issues along this avenue. 
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